Análisis de Componentes Principales

Asignatura

Análisis Multivariante

Profesor: Laura Durán

Alumnos: Steven Allus, Antonio López


1. Descripción de la actividad

En esta actividad, exploraremos el Análisis de Componentes Principales (PCA) aplicando técnicas estadísticas en R.

La tarea consiste en llevar a cabo un Análisis de Componentes Principales a partir de un estudio en el que se ha medido la expresión de los genes sobre distintos tumores. Los datos están almacenados en el archivo genes_tumores.csv que se puede descargar desde el aula virtual. Esta tabla contiene los valores de expresión de 12 genes medidos en 10 tumores distintos.

El objetivo del análisis será averiguar si de los 10 tumores analizados existen algunos con características comunes en cuanto a la expresión de sus genes. Los resultados del análisis deberán representarse en un biplot en el que las variables (genes) se representen como vectores y los individuos (tumores) como puntos.

Pueden seguirse los pasos explicados en clase u otros alternativos que conduzcan a las mismas conclusiones.

El trabajo es en parejas

Aspectos formales de las tareas

Hay que entregar un informe con el código en R utilizado, las explicaciones pertinentes para cada paso del análisis y las conclusiones derivadas del análisis. Se recomienda usar RMarkdown.


2. Definición Componentes Principales

Los Componentes Principales son combinaciones lineales de las variables originales diseñadas para maximizar la varianza explicada.

Datos con los que trabajar: Diez tumores con datos de doce genes.

Tenemos que reducir las variables de los tumores para representarlas visualmente siempre perdiendo información.

A partir de la combinación lineal de las doce variables obtendremos una reducción de estas que nos permitirá representarlas gráficamente.

Estas variables representadas gráficamente se llaman Componentes Principales.


3. Ejemplo representación tres variables

Vamos a representar un eje de coordenadas utilizando subconjuntos de datos de los genes. Cada conjunto incluye tres genes, lo que permite una visualización tridimensional intuitiva. Los datos están representados con un gráfico interactivo que muestra cómo se distribuyen las muestras de tumores en función de las combinaciones lineales de los genes. Esta representación facilita observar similitudes o diferencias entre tumores en relación con los genes seleccionados.


4. Código R

4.1 Instalación librerías, carga de datos y muestra de datos

# Establecer el repositorio de CRAN para evitar errores
options(repos = c(CRAN = "https://cloud.r-project.org"))

# Seleccionar directorio de trabajo
setwd("/Users/stevenallus/Downloads")

# Instalar y cargar paquetes necesarios
if (!require(boot)) install.packages("boot")
library(boot)
paquetes <- c("logisticPCA", "FactoMineR", "factoextra", "plotly", "readxl")
paquetes_faltantes <- paquetes[!paquetes %in% installed.packages()[,"Package"]]
if (length(paquetes_faltantes) > 0) install.packages(paquetes_faltantes)

# Cargar las librerías
lapply(paquetes, library, character.only = TRUE)
## [[1]]
## [1] "logisticPCA" "boot"        "stats"       "graphics"    "grDevices"  
## [6] "utils"       "datasets"    "methods"     "base"       
## 
## [[2]]
##  [1] "FactoMineR"  "logisticPCA" "boot"        "stats"       "graphics"   
##  [6] "grDevices"   "utils"       "datasets"    "methods"     "base"       
## 
## [[3]]
##  [1] "factoextra"  "ggplot2"     "FactoMineR"  "logisticPCA" "boot"       
##  [6] "stats"       "graphics"    "grDevices"   "utils"       "datasets"   
## [11] "methods"     "base"       
## 
## [[4]]
##  [1] "plotly"      "factoextra"  "ggplot2"     "FactoMineR"  "logisticPCA"
##  [6] "boot"        "stats"       "graphics"    "grDevices"   "utils"      
## [11] "datasets"    "methods"     "base"       
## 
## [[5]]
##  [1] "readxl"      "plotly"      "factoextra"  "ggplot2"     "FactoMineR" 
##  [6] "logisticPCA" "boot"        "stats"       "graphics"    "grDevices"  
## [11] "utils"       "datasets"    "methods"     "base"
# Importar datos desde un archivo CSV
genes_tumores <- read.csv("genes_tumores.csv")
genes_tumores <- genes_tumores[, -1] # Excluir la primera columna si no es relevante
head(genes_tumores)
##        gene1      gene2      gene3      gene4      gene5       gene6
## 1  0.9323190  1.2113278  0.5383686 -0.3804224 -0.5994117 -0.46955658
## 2  1.0613959  0.7762260  0.5098463  0.1574786 -0.4035551 -0.50984955
## 3  0.4637662  1.1793956  1.1355446  0.1548933 -0.7319309 -0.42653172
## 4 -1.6228992 -0.7777483 -0.7868679 -1.1771294 -0.3007739  0.11591201
## 5 -1.1035586 -0.9511971 -0.9131699 -1.0484974 -0.3974149 -0.35492610
## 6 -1.0950318 -0.9722049 -0.9504966 -0.8261171 -0.5278552 -0.06271155
##         gene7       gene8      gene9    gene10    gene11    gene12
## 1 -1.06758821 -1.23610151  0.1145246 0.8975893 1.4624217 1.1259876
## 2 -1.31883371 -0.87248832 -0.2017863 1.5210825 1.1883528 0.9711843
## 3 -1.03071756 -1.05634140 -0.2747214 1.0444794 0.7132720 0.7550505
## 4  0.19523024  0.15654077  0.7884622 0.9194994 1.0616514 0.9178455
## 5  0.23230945 -0.22601482  1.0759271 1.1158295 0.8121759 1.1303828
## 6 -0.04847066  0.05124849  1.3769491 0.8643747 0.9175155 0.6019269

4.2 Escalado variables, reducción de dimensiones

Este código escala todas las variables para conseguir las componentes principales y nos muestra un gráfico.

El escalado (normalización) asegura que todas las variables tengan la misma importancia en el análisis, independientemente de sus unidades o magnitudes. Esto se hace transformando cada variable para que tenga una media de 0 y una desviación estándar de 1.

Sin el escalado, las variables con valores más grandes (como pesos en kilogramos) dominarían a las variables con valores más pequeños (como proporciones o porcentajes).

Las componentes principales son combinaciones lineales de las variables originales, y el escalado garantiza que estas combinaciones no estén sesgadas por las escalas originales de las variables. Esto permite que el análisis se enfoque únicamente en las relaciones entre las variables y no en sus magnitudes.

# PCA -> datos escalados, con plot:
      
      pca_genes <- PCA(genes_tumores, scale.unit = T, graph = T)

Gráfico de PCA mostrando la reducción de dimensiones.

Conclusiones:

  • Dim 1 (PC1): Explica el 68.19% de la varianza.
  • Dim 2 (PC2): Explica el 28.26% de la varianza.
  • Componentes principales = 96.45% de la varianza total.

Esto significa que los dos primeros componentes principales explican el 96.45% de la información contenida en los datos originales.


4.3 Eigenvalues

Ordena las dimensiones de mayor a menor valor de eigenvalue.

Eigenvalues o valores propios: mide la cantidad de varianza explicada por esa componente.

# str(pca_coronaria) para ver la estructura del objeto
    summary(pca_genes)
## 
## Call:
## PCA(X = genes_tumores, scale.unit = T, graph = T) 
## 
## 
## Eigenvalues
##                        Dim.1   Dim.2   Dim.3   Dim.4   Dim.5   Dim.6   Dim.7
## Variance               8.183   3.391   0.132   0.095   0.070   0.053   0.050
## % of var.             68.191  28.260   1.098   0.789   0.579   0.443   0.417
## Cumulative % of var.  68.191  96.451  97.549  98.338  98.918  99.361  99.778
##                        Dim.8   Dim.9
## Variance               0.017   0.009
## % of var.              0.144   0.079
## Cumulative % of var.  99.921 100.000
## 
## Individuals
##            Dist    Dim.1    ctr   cos2    Dim.2    ctr   cos2    Dim.3    ctr
## 1      |  2.713 | -1.819  4.043  0.450 |  1.903 10.682  0.492 | -0.141  1.513
## 2      |  2.594 | -1.543  2.910  0.354 |  1.928 10.962  0.552 |  0.485 17.858
## 3      |  2.494 | -1.071  1.401  0.184 |  2.085 12.822  0.699 | -0.654 32.420
## 4      |  3.363 | -2.178  5.798  0.419 | -2.485 18.207  0.546 | -0.220  3.678
## 5      |  3.464 | -2.656  8.622  0.588 | -2.160 13.763  0.389 |  0.186  2.614
## 6      |  3.391 | -2.363  6.824  0.486 | -2.352 16.316  0.481 |  0.096  0.694
## 7      |  4.218 |  4.100 20.541  0.945 | -0.717  1.518  0.029 | -0.470 16.784
## 8      |  4.199 |  4.168 21.230  0.985 | -0.271  0.216  0.004 | -0.032  0.079
## 9      |  4.711 |  4.662 26.566  0.980 | -0.215  0.136  0.002 |  0.515 20.123
## 10     |  2.697 | -1.300  2.065  0.232 |  2.284 15.377  0.717 |  0.236  4.238
##          cos2  
## 1       0.003 |
## 2       0.035 |
## 3       0.069 |
## 4       0.004 |
## 5       0.003 |
## 6       0.001 |
## 7       0.012 |
## 8       0.000 |
## 9       0.012 |
## 10      0.008 |
## 
## Variables (the 10 first)
##           Dim.1    ctr   cos2    Dim.2    ctr   cos2    Dim.3    ctr   cos2  
## gene1  |  0.530  3.434  0.281 |  0.808 19.260  0.653 |  0.198 29.684  0.039 |
## gene2  |  0.577  4.066  0.333 |  0.800 18.876  0.640 | -0.093  6.571  0.009 |
## gene3  |  0.595  4.321  0.354 |  0.787 18.268  0.620 | -0.100  7.620  0.010 |
## gene4  |  0.910 10.117  0.828 |  0.372  4.079  0.138 |  0.019  0.279  0.000 |
## gene5  |  0.949 11.014  0.901 | -0.220  1.433  0.049 |  0.191 27.663  0.036 |
## gene6  |  0.917 10.278  0.841 | -0.378  4.217  0.143 | -0.040  1.211  0.002 |
## gene7  |  0.707  6.114  0.500 | -0.674 13.390  0.454 | -0.102  7.833  0.010 |
## gene8  |  0.781  7.449  0.610 | -0.611 11.002  0.373 |  0.092  6.401  0.008 |
## gene9  | -0.867  9.185  0.752 | -0.449  5.952  0.202 |  0.051  1.980  0.003 |
## gene10 | -0.963 11.328  0.927 |  0.198  1.154  0.039 |  0.064  3.090  0.004 |

Los eigenvalues miden la varianza explicada por cada componente.


4.4 Componentes principales Dim1 y Dim2

plot(
    pca_genes$ind$coord[, 1], # Coordenadas en Dim1
      pca_genes$ind$coord[, 2], # Coordenadas en Dim2
      xlab = paste0("Dim 1 (", round(pca_genes$eig[1, 2], 2), "% varianza)"), 
      ylab = paste0("Dim 2 (", round(pca_genes$eig[2, 2], 2), "% varianza)"), 
      main = "Proyección en Dim1 y Dim2", # Título del gráfico
      pch = 19, # Tipo de punto
      col = "blue", # Color de los puntos
      xlim = c(min(pca_genes$ind$coord[, 1]) - 1, max(pca_genes$ind$coord[, 1]) + 1), 
      ylim = c(min(pca_genes$ind$coord[, 2]) - 1, max(pca_genes$ind$coord[, 2]) + 1)  )
    
    # Eje coordenadas
    abline(h = 0, v = 0, col = "black", lty = 2) 
    
    # Etiquetas 
    text(
      pca_genes$ind$coord[, 1],
      pca_genes$ind$coord[, 2],
      labels = rownames(pca_genes$ind$coord),
      pos = 4, # Posición del texto
      cex = 0.7 # Tamaño del texto
    )

Resultados —

Proyección de los datos en Dim1 y Dim2.

4.5 Varianza de cada dimensión

# Varianza explicada por cada dimensión:
fviz_eig(pca_genes, addlabels = TRUE, ylim = c(0, 100))


4.6 Varianza acumulada en cada dimensión

library(ggplot2)
library(plotly)

# Crear el data frame
varianza <- data.frame(
  Dimensión = 1:length(pca_genes$eig[, 1]),
  Varianza = pca_genes$eig[, 2],
  Acumulada = cumsum(pca_genes$eig[, 2])
)

# Gráfico con ggplot2
graf_varianza <- ggplot(varianza, aes(x = Dimensión)) +
  geom_bar(aes(y = Varianza), stat = "identity", fill = "steelblue") +
  geom_line(aes(y = Acumulada), color = "red", linewidth = 1) +
  geom_point(aes(y = Acumulada), color = "red", size = 2) +
  labs(
    title = "Varianza Explicada y Acumulada",
    x = "Dimensión",
    y = "Porcentaje de Varianza"
  ) +
  theme_minimal()

# Exportar el gráfico a plotly
ggplotly(graf_varianza)
## Warning: 'bar' objects don't have these attributes: 'mode'
## Valid attributes include:
## '_deprecated', 'alignmentgroup', 'base', 'basesrc', 'cliponaxis', 'constraintext', 'customdata', 'customdatasrc', 'dx', 'dy', 'error_x', 'error_y', 'hoverinfo', 'hoverinfosrc', 'hoverlabel', 'hovertemplate', 'hovertemplatesrc', 'hovertext', 'hovertextsrc', 'ids', 'idssrc', 'insidetextanchor', 'insidetextfont', 'legendgroup', 'legendgrouptitle', 'legendrank', 'marker', 'meta', 'metasrc', 'name', 'offset', 'offsetgroup', 'offsetsrc', 'opacity', 'orientation', 'outsidetextfont', 'selected', 'selectedpoints', 'showlegend', 'stream', 'text', 'textangle', 'textfont', 'textposition', 'textpositionsrc', 'textsrc', 'texttemplate', 'texttemplatesrc', 'transforms', 'type', 'uid', 'uirevision', 'unselected', 'visible', 'width', 'widthsrc', 'x', 'x0', 'xaxis', 'xcalendar', 'xhoverformat', 'xperiod', 'xperiod0', 'xperiodalignment', 'xsrc', 'y', 'y0', 'yaxis', 'ycalendar', 'yhoverformat', 'yperiod', 'yperiod0', 'yperiodalignment', 'ysrc', 'key', 'set', 'frame', 'transforms', '_isNestedKey', '_isSimpleKey', '_isGraticule', '_bbox'

4.7 BIPLOT simultáneo variables e individuos

# Librerías necesarias
library(FactoMineR)
library(factoextra)
library(plotly)

# Crear el biplot con FactoMineR
biplot_pca <- fviz_pca_biplot(
  pca_genes,
  repel = TRUE,            # Evitar superposición de etiquetas
  col.var = "blue",        # Color de las flechas de variables
  col.ind = "black",       # Color de los puntos de individuos
  arrowsize = 1,           # Tamaño de las flechas
  pointsize = 3,           # Tamaño de los puntos
  geom = c("point", "text"),
  title = "Biplot: Genes (Flechas) y Tumores (Puntos)"
)

# Exportar el biplot a plotly para interactividad
ggplotly(biplot_pca)
## Warning in geom2trace.default(dots[[1L]][[1L]], dots[[2L]][[1L]], dots[[3L]][[1L]]): geom_GeomTextRepel() has yet to be implemented in plotly.
##   If you'd like to see this geom implemented,
##   Please open an issue with your example code at
##   https://github.com/ropensci/plotly/issues
## Warning in geom2trace.default(dots[[1L]][[1L]], dots[[2L]][[1L]], dots[[3L]][[1L]]): geom_GeomTextRepel() has yet to be implemented in plotly.
##   If you'd like to see this geom implemented,
##   Please open an issue with your example code at
##   https://github.com/ropensci/plotly/issues

Conclusiones

  1. Interpretación de los genes en las componentes principales:
    • En Dim1 (PC1), los genes con mayor contribución son:
      • gene10 (-0.963): Este gen tiene una contribución significativa y opuesta a otros genes en el eje PC1, lo que indica que es clave para diferenciar tumores con características genéticas opuestas.
      • gene5 (0.949): Su alta contribución positiva sugiere una fuerte relación con el eje PC1, lo que lo convierte en un gen fundamental para discriminar tumores en esta dimensión.
    • En Dim2 (PC2), los genes destacados son:
      • gene1 (0.808): Este gen tiene un peso importante en PC2, influyendo en la diferenciación de tumores que comparten ciertas similitudes genéticas.
      • gene3 (0.787): Similar al gene1, este gen juega un papel relevante para identificar subgrupos de tumores, probablemente aquellos con una regulación diferencial específica.
  2. Análisis de los tumores en el biplot:
    • Tumores agrupados:
      • Los tumores G (PC1=4.100, PC2=0.717), H (PC1=4.168, PC2=0.271) y I (PC1=4.662, PC2=0.214) tienen posiciones cercanas en el espacio de las componentes principales. Esto indica similitudes significativas en la expresión genética, lo que sugiere que estos tumores podrían compartir patrones de comportamiento o ser tratados bajo estrategias similares en investigaciones futuras.
    • Tumores atípicos:
      • Los tumores D (PC1=-2.178, PC2=2.485) y E (PC1=-2.656, PC2=2.160) se posicionan lejos de los demás, lo que sugiere perfiles genéticos únicos o características que podrían requerir atención especializada en futuros análisis.
  3. Implicaciones:
    • Estos resultados permiten identificar tumores con patrones genéticos similares, lo cual podría ser de gran utilidad en la clasificación y diseño de tratamientos personalizados. Por ejemplo, tumores como G, H e I podrían investigarse como un subgrupo homogéneo, mientras que los tumores atípicos D y E requerirían un enfoque más específico.

4.9 Extra Metodo de Validación

Vamos a validar la estabilidad y consistencia del análisis realizado. Para ello haremos empleo la técnica de remuestreo (bootstrap) para evaluar la estabilidad de los componentes principales. Aquí está la propuesta:

# Función para remuestrear los datos y recalcular el PCA
pca_bootstrap <- function(data, indices) {
  # Subconjunto de datos remuestreados
  data_boot <- data[indices, ]
  
  # Realizar PCA con los datos remuestreados
  pca_result <- PCA(data_boot, scale.unit = TRUE, graph = FALSE)
  
  # Concatenar las coordenadas de PC1 y PC2 en un vector
  c(pca_result$ind$coord[, 1], pca_result$ind$coord[, 2])
}

# Configurar remuestreo bootstrap
set.seed(123) # Para reproducibilidad
bootstrap_results <- boot(
  data = genes_tumores,           # Datos originales
  statistic = pca_bootstrap,      # Función definida arriba
  R = 1000                        # Número de remuestreos
)

# Resumen de resultados
print(bootstrap_results)
## 
## ORDINARY NONPARAMETRIC BOOTSTRAP
## 
## 
## Call:
## boot(data = genes_tumores, statistic = pca_bootstrap, R = 1000)
## 
## 
## Bootstrap Statistics :
##        original     bias    std. error
## t1*  -1.8188739  2.0286898    2.995240
## t2*  -1.5431329  1.6136273    2.914554
## t3*  -1.0708431  0.8549935    2.761629
## t4*  -2.1782566  2.0206275    2.808105
## t5*  -2.6561326  2.7046518    2.950829
## t6*  -2.3630028  2.3056353    2.877474
## t7*   4.0998488 -4.1402529    2.853162
## t8*   4.1679928 -4.0180110    2.850127
## t9*   4.6624718 -4.7633048    2.823456
## t10* -1.3000715  1.3933435    2.936113
## t11*  1.9032783 -1.8512380    1.749042
## t12*  1.9280936 -1.9604068    1.804571
## t13*  2.0852384 -2.1047527    1.862594
## t14* -2.4848021  2.5909666    1.831133
## t15* -2.1603981  2.1103791    1.762714
## t16* -2.3522557  2.3860561    1.771807
## t17* -0.7174422  0.7045590    1.848779
## t18* -0.2705048  0.2949778    1.787310
## t19* -0.2147895  0.1865351    1.854436
## t20*  2.2835821 -2.3570762    1.856363
# Visualización de la variación en las componentes principales
plot(bootstrap_results, index = 1, main = "Variación en PC1 (Bootstrap)")

plot(bootstrap_results, index = 2, main = "Variación en PC2 (Bootstrap)")


De los resultados del análisis de bootstrap para los tumores (t1* a t19*), se pueden extraer las siguientes conclusiones basadas en los valores presentados y los gráficos (histograma y QQ-plot):


Conclusiones:

1. Estabilidad de las Coordenadas (Bias y Std. Error):

  • Bias: Representa el promedio de las diferencias entre las coordenadas originales y las obtenidas en los remuestreos.
    • Tumores como t7, t8 y t9* tienen un bias alto (negativo), lo que indica una mayor sensibilidad en sus coordenadas respecto al remuestreo. Esto podría deberse a su posición más extrema en el espacio de las componentes principales.
    • Tumores como t18* y t19* presentan bias más cercanos a 0, lo que sugiere una mayor estabilidad en sus posiciones.
  • Std. Error: Representa la variabilidad en las coordenadas a través de los remuestreos.
    • Tumores como t11* a t19* tienen un menor error estándar (alrededor de 1.7 - 1.8), lo que indica que sus posiciones son más consistentes.
    • Tumores como t1* a t10* presentan mayor error estándar (alrededor de 2.8 - 3.0), indicando más variación en los remuestreos y menor estabilidad en sus posiciones.

2. Distribución de las Coordenadas (Gráficos):

  • Histograma:
    • La distribución de las coordenadas generadas por bootstrap es aproximadamente simétrica, aunque hay casos con colas largas, lo que indica posible sensibilidad a valores extremos.
  • QQ-plot:
    • La mayoría de los puntos siguen la línea de normalidad, pero los extremos se desvían significativamente. Esto sugiere que algunos tumores tienen posiciones menos robustas (posiblemente los outliers identificados en el análisis).

3. Interpretación de Tumores Específicos:

  • Tumores t7, t8 y t9*:
    • Presentan bias negativos altos y altos errores estándar, lo que sugiere que son los más sensibles al remuestreo. Esto está alineado con su posición atípica en el biplot original.
  • Tumores t11* a t19*:
    • Menor bias y std. error indican que son más consistentes en sus posiciones y menos sensibles al remuestreo, fortaleciendo su agrupación observada en el biplot.

Conclusión Final:

El análisis bootstrap confirma que: 1. Los patrones de agrupación observados en el biplot son consistentes, especialmente para tumores con menor bias y std. error. 2. Algunos tumores atípicos (como t7, t8, t9* y t14*), aunque claramente identificados en el biplot, son más sensibles al remuestreo, lo que sugiere que sus posiciones deben interpretarse con precaución. 3. Tumores centrales en el biplot tienen posiciones robustas y estables, reforzando su relevancia en el análisis.